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Abstract. In this paper we study the formal algebraic structure underlying 
the intrinsic classification algorithm, recently introduced by Hadani, Shkol- 
nisky, Singer and Zhao, for classifying noisy projection images of similar view- 
ing directions in three-dimensional cryo-electron microscopy (cryo-EM). This 
preliminary classification is of fundamental importance in determining the 
three-dimensional structure of macromolecules from cryo-EM images. Inspect- 
ing this algebraic structure we obtain a conceptual explanation for the admis- 
sibility (correctness) of the algorithm and a proof of its numerical stability, 
thus putting it on firm mathematical grounds. The proof relies on studying 
the spectral properties of an integral operator of geometric origin on the two- 
dimensional sphere, called the localized parallel transport operator. Along 
the way, we continue to develop the representation theoretic setup for three- 
dimensional cryo-EM that was initiated in [12) . 



0. Introduction 

The goal in cryo-EM is to determine the three-dimensional structure of a mol- 
ecule from noisy projection images taken at unknown random orientations by an 
electron microscope, i.e., a random Computational Tomography (CT). Determining 
three-dimensional structures of large biological molecules remains vitally important, 
as witnessed, for example, by the 2003 Chemistry Nobel Prize, co-awarded to R. 
MacKinnon for resolving the three-dimensional structure of the Shaker K-f channel 
protein [H [13], and by the 2009 Chemistry Nobel Prize, awarded to V. Ramakr- 
ishnan, T. Steitz and A. Yonath for studies of the structure and function of the 
ribosome. The standard procedure for structure determination of large molecules 
is X-ray crystallography. The challenge in this method is often more in the crys- 
tallization itself than in the interpretation of the X-ray results, since many large 
molecules, including various types of proteins have so far withstood all attempts to 
crystallize them. 

Cryo-EM is an alternative approach to X-ray crystallography. In this approach, 
samples of identical molecules are rapidly immobilized in a thin layer of vitreous 
ice (this is an ice without crystals) . The cryo-EM imaging process produces a large 
collection of tomographic projections, corresponding to many copies of the same 
molecule, each immobilized in a different (yet unknown) orientation. The intensity 
of the pixels in a given projection image is correlated with the line integrals of the 
electric potential induced by the molecule along the path of the imaging electrons 
(see Figure [T]) . The goal is to reconstruct the three-dimensional structure of the 
molecule from such a collection of projection images. The main problem is that the 
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highly intense electron beam damages the molecule and, therefore, it is problematic 
to take projection images of the same molecule at known different directions as in 
the case of classical CtQ. In other words, a single molecule is imaged only once, 
rendering an extremely low signal-to- noise ratio (SNR), mostly due to shot noise 
induced by the maximal allowed electron dose. 




Figure 1. Schematic drawing of the imaging process: every pro- 
jection image corresponds to some unknown spatial orientation of 
the molecule. 



0.1. Mathematical model. Instead of thinking of a multitude of molecules immo- 
bilized in various orientations and observed by a microscope held in a fixed position, 
it is more convenient to think of a single molecule, observed by an electron micro- 
scope from various orientations. Thus, an orientation describes a configuration of 
the microscope instead of that of the molecule. 

Let {V, {■,■)) be an oriented three-dimensional Euclidean vector space. The 
reader can take V to be R'' and (•, •) to be the standard inner product. Let 
X — Ft (V) be the oriented frame manifold associated to V; a point x E X is 
an orthonormal basis x = (61,62,63) of V compatible with the orientation. The 
third vector 63 is distinguished, denoted by tt (x) and called the viewing direction. 
More concretely, if we identify V with R'^, then a point in X can be thought of as 
a matrix belonging to the special orthogonal group SO (3), whose first, second and 
third columns are the vectors 61,62 and 63 respectively. 

Using this terminology, the physics of cryo-EM is modeled as follows: 

• The molecule is modeled by a real valued function cj) : V ^ M., describing 
the electromagnetic potential induced from the charges in the molecule. 



We remark that there are other methods like single-or multi-ajcis tilt EM tomogrophy, where 
several lower dose/higher noise images of a single molecule are taken from known directions. These 
methods are used for example when one has an organic object in vitro or a collection of different 
objects in the sample. There is a rich literature for this field starting with the work of Crowther, 
DeRosier and Klug in the early 1960s. 
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• A spatial orientation of the microscope is modeled by an ortlionormal frame 
X dz X . The third vector tt (x) is the viewing direction of the microscope 
and the plane spanned by the first two vectors ei and 62 is the plane of the 
camera equipped with the coordinate system of the camera (see Figure [2]) • 

• The projection image obtained by the microscope, when observing the mol- 
ecule from a spatial orientation a; is a real valued function / : — >■ M, given 
by the X-ray projection along the viewing direction: 

I {P, q) = Xray^(^)0 {p, q) ^ J cj) (pei + qe2 + teg) dr. 

teR 

for every (p, g) e M'^ . 

The data collected from the experiment is a set consisting of A'' projection images 
V = {/i, ..,In}- Assuming that the potential function cj) is generico, in the sense 
that, each image li € V can originate from a unique frame Xi € X , the main problem 
of cryo-EM is to reconstruct the (unique) unknown frame Xi € X associated with 
each projection image li S V. 




Figure 2. A frame x — (ei, 62, 63) modeling the orientation of the 
electron microscope, where tt [x) — 63 is the viewing direction and 
the pair (ei, 62) establishes the coordinates of the camera. 

0.2. Class averaging. As projection images in cryo-EM have extremely low SNR0 
(see Figure[3]), a crucial initial step in all reconstruction methods is "class averaging" 
[8]. Class averaging is the grouping of a large data set of noisy raw projection 
images into clusters, such that images within a single cluster have similar viewing 
directions. Averaging rotationally aligned noisy images within each cluster results 
in "class averages"; these are images that enjoy a higher SNR and are used in 
later cryo-EM procedures such as the angular reconstitution procedure, IQ^, that 



This assumption about the potential can be omitted in the context of the class averaging 
algorithm presented in this paper. In particular, the algorithm can be applied to potentials 
describing molecules with symmetries which do not satisfy the "generic" assumption. 

■^SNR stands for Signal to Noise Ratio, which is the ratio between the squared norm of the 
signal and the squared norm of the noise. 
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requires better quality images. Finding consistent class averages is challenging due 
to the high level of noise in the raw images. 




Figure 3. The left most image is a clean simulated projection 
image of the E.coli SOS ribosomal subunit. The other three images 
are real electron microscope images of the same subunit. 



The starting point for the classification is the idea that visual similarity between 
projection images suggests vicinity between viewing directions of the corresponding 
(unknown) frames. The similarity between images li and Ij is measured by their 
invariant distance (introduced in [111) which is the Euclidean distance between the 
images when they are optimally aligned with respect to in-plane rotations, namely 

(0.1) (/.,/,)= min ||i?(g)/,-/, ||, 

geSO(2) 

where 

for any function / : — >• R. 

One can choose some threshold value e, such that d < e is indicative that 

perhaps the corresponding frames Xi and xj have nearby viewing directions. The 
threshold e defines an undirected graph G — (Vertices, Edges) with vertices labeled 
by numbers 1,..,A'^ and an edge connecting vertex i with vertex j if and only if 
the invariant distance between the projection images Jj and Ij is smaller then e, 
namely 

{i,j} e Edges d < e. 

In an ideal noiseless world, the graph G acquires the geometry of the unit sphere 
S (V), namely, two images are connected by an edge if and only if their correspond- 
ing viewing directions are close on the sphere, in the sense that they belong to some 
small spherical cap of opening angle a = a (e) . 

However, the real world is far from ideal as it is governed by noise; hence, it 
often happens that two images of completely different viewing directions have small 
invariant distance. This can happen when the realizations of the noise in the two 
images match well for some random in-plane rotation, leading to spurious neighbor 
identification. Therefore, the nai've approach of averaging the rotationally aligned 
nearest neighbor images can sometimes yield a poor estimate of the true signal in 
the reference image. 

To summarize: From this point of view, the main problem is to distinguish 
the good edges from the bad ones in the graph G, or, in other words, to distin- 
guish the true neighbors from the false ones (called outliers). The existence of 



THE CLASS AVERAGING PROBLEM 



5 



outliers is the reason why the classification problem is non-trivial. We emphasize 
that without excluding the outliers, averaging rotationally aligned images of small 
invariant distance (j0.1[) yields poor estimate of the true signal, rendering the prob- 
lem of three-dimensional reconstruction from cryo-EM images non-feasible. In this 
respect, the class averaging problem is of fundamental importance. 

0.3. Main results. In 2 , we introduced a novel algorithm, referred to in this 
paper as the intrinsic classification algorithm, for classifying noisy projection images 
of similar viewing directions. The main appealing property of this new algorithm is 
its extreme robustness to noise and to presence of outliers; in addition, it also enjoys 
efficient time and space complexity. These properties are explained thoroughly in 
[2], which includes also a large number of numerical experiments. 

In this paper we study the formal algebraic structure that underlies the intrinsic 
classification algorithm. Inspecting this algebraic structure we obtain a conceptual 
explanation for the admissibility (correctness) of the algorithm and a proof of its 
numerical stability, thus putting it on firm mathematical grounds. The proof relies 
on the study of a certain integral operator Th on X , of geometric origin, called the 
localized parallel transport operator. Specifically: 

• Admissibility amounts to the fact that the maximal eigenspace of Th is 
a three-dimensional complex Hermitian vector space and that there is a 
canonical identification of Hermitian vector spaces between this eigenspace 
and the complexified vector space W = CV. 

• Numerical stability amounts to the existence of a spectral gap which sep- 
arates the maximal eigenvalue of from the rest of the spectrum, which 
enables one to obtain a stable numerical approximation of the correspond- 
ing maximal eigenspace and of other related geometric structures. 

The main technical result of this paper is a complete description of the spectral 
properties of the localized parallel transport operator. Along the way, we continue 
to develop the mathematical set-up for cryo-EM that was initiated in [12], thus 
further elucidating the central role played by representation theoretic principles in 
this scientific discipline. 

The remainder of the introduction is devoted to a detailed description of the 
intrinsic classification algorithm and to an explanation of the main ideas and results 
of this paper. 

0.4. Transport data. A preliminary step is to extract certain geometric data from 
the set of projection images, called (local) empirical transport data. 

When computing the invariant distance between images li and Ij we also record 
the rotation matrix in 5*0 (2) that realizes the minimum in (10. ip and denote this 
special rotation by T (i, j), that is 

(0.2) r(z,j) = argmin||i?(g)/,-/,||. 

seS'0(2) 

noting that, 

(0.3) f{J,^)^f{^,J)-\ 

The main observation is that in an ideal noiseless world the rotation T {i,j) can 
be interpreted as a geometric relation between the corresponding frames Xi and Xj , 
provided the invariant distance between the corresponding images is small. This 
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relation is expressed in terms of parallel transport on the sphere, as follows: define 
the rotation 

l,sin(%) cos(%);' 
as the unique solution of the equation 

(0.4) Xi <] T {xi, Xj) — ,n.(j.j)a;j , 

where i7r(3;i).7r(xj) is the parallel transport along the unique geodesic on the sphere 
connecting the points (xj) with n (xi) or, in other words, it is the rotation in 
SO {V) that takes the vector tt (xj) to tt {xi) along the shortest path on the sphere 
and the action < is defined by 

(si'nl^j cos(^)) = (cos (61) ei + sin (6*) 62, -sin (6*) ei + cos (6*) 62, 63), 

for every x = (61, 62, 63). The precise statement is that the rotation T approx- 
imates the rotation T {xi,Xj) when £ Edges. This geometric interpretation of 
the rotation T is suggested from a combination of mathematical and empirical 
considerations that we proceed to explain. 

• On the mathematical side: the rotation T{xi,Xj) is the unique rotation of 
the frame Xi around its viewing direction tt {xi), minimizing the distance to 
the frame Xj. This is a standard fact from differential geometry (a direct 
proof of this statement appears in '^). 

• On the empirical side: if the function (p is "nice", then the optimal align- 
ment T{i,j) of the projection images is correlated with the optimal align- 
ment T{xi,Xj) of the corresponding frames. This correlation of course 
improves as the distance between tt (xi) and tt (xj) becomes smaller. A 
quantitative study of the relation between T(i,j) and T(xi,Xj) involves 
considerations from image processing thus it is beyond the scope of this 
paper. 

To conclude, the "empirical" rotation T{i,j) approximates the "geometric" ro- 
tation T [xi, Xj) only when the viewing directions tt (xi) and tt (xj) are close, in the 
sense that they belong to some small spherical cap of opening angle a. The later 
"geometric" condition is correlated with the "empirical" condition that the corre- 
sponding images li and Ij have small invariant distance. When tt (xi) and tt (xj) 
are far from each other, the rotation T{i,j) is not related any longer to parallel 
transportation on the sphere. For this reason, we consider only rotations T (j, j) for 
which {i,j} S Edges and call this collection the (local) empirical transport data. 

0.5. The intrinsic classification algorithm. The intrinsic classification algo- 
rithm accepts as an input the empirical transport data {T {i,j) : £ Edges} and 
produces as an output the Euclidean inner products {(tt (xi) , tt {xj)) : i,j = 1, ..N}. 
Using these inner products, one can identify the true neighbors in the graph G, as 
the pairs {i,j} £ Edges for which the inner product {tt {xi) ,tt {xj j) is close to 1. 
The formal justification of the algorithm requires the empirical assumption that the 
frames Xi, i = 1, ..,iV are uniformly distributed in the frame manifold X, accord- 
ing to the unique normalized Haar measure on X. This assumption corresponds 
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to the situation where the orientations of the molecules in the ice are distributed 
independently and uniformly at random. 

The main idea of the algorithm is to construct an intrinsic model, denoted by 
Wat, of the Hermitian vector space W = CV which is expressed solely in terms of 
the empirical transport data. 

The algorithm proceeds as follows: 

Stepl (Ambient Hilbert space): consider the standard iV-dimensional Hilbert 
space 

Hn = C^. 

Step 2 (Self adjoint operator): identify with C and consider each rotation 
T{i,j) as a complex number of unit norm. Define the N x N complex matrix 

by putting the rotation T{i,j) in the entry. Notice that the matrix is 

self-adjoint by ()0.3|) . 

Step 3 (Intrinsic model): the matrix Tn induces a spectral decomposition 

HAf=0HjV (A). 
A 

Theorem 1. There exists a threshold Xq such that 

dim ^^(A) = 3. 

A>Ao 

Define the Hermitian vector space 

Wat = (A) . 

A>Ao 

Step 4 (Computation of the EucUdean inner products): the Euclidean 
inner products {(tt (xi) , tt {xj)) : i,j = 1, ..N} are computed from the vector space 
Wiv, as follows: for every i = 1,..,N, denote by G Wat the vector 

-y273-pr*(l), 

where pr^ : Wjv — > C is the projection on the ith component and pr* : C — > Wjv 
is the adjoint map. In addition, for every frame x € X , x = (61,62,63), denote by 
Sx & W the (complex) vector 61 — ie2- 

The upshot is that the intrinsic vector space W^r consisting of the collection of 
vectors ip^ g W^r, i — I, .., N is (approximateljQ) isomorphic to the extrinsic vector 
space W consisting of the collection of vectors S^^ & W , i — 1, .., N , where Xi is the 
frame corresponding to the image 7^, for every i — 1,..,N. This statement is the 
content of the following theorem: 

Theorem 2. There exists a unique (approximated) isomorphism tn '■ W ^ Wat 
of Hermitian vector spaces such that 

TN {Sx,) = (fi, 

for every i = 1, .., N . 



This approximation improves as grows. 
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The above theorem enables us to express, in intrinsic terms, the Euchdean inner 
products between the viewing directions, as follows: starting with the following 
identity from linear algebra (that will be proved in the sequel): 

(0.5) {^{x),ii{y))^\{5,,6y)\^l, 

for every pair of frames x,y € X, where (•, •) is the Euclidean product on V and 
(•, •) is the Hermitian product onW = CV , induced from (•, •), given by 

{u + iv, u' + iv') ~ {u, v) + {v, v') — i {u, v') + i {v, u') , 

we obtain the following relation: 

(0.6) {-Kixi) ,tt{xj)) = \ {ip^,Lp^) \ ~ 1, 

for every i, j — 1, .., N. In the derivation of Relation (|0.6p from Relation (|0.5p we 
use Theorem [2] Notice that Relation (|0.6[) implies that although we do not know 
the frame associated with every projection image, we still are able to compute the 
inner product between every pair of such frames from the intrinsic vector space Wjv 
which, in turns, can be computed from the images. 

0.6. Structure of the paper. The paper consists of three sections besides the 
introduction. 

• In Section [l] we begin by introducing the basic analytic setup which is 
relevant for the class averaging problem in cryo-EM. Then, we proceed to 
formulate the main results of this paper, which are: a complete description 
of the spectral properties of the localized parallel transport operator (The- 
orem [3]) , the spectral gap property (Theorem 3]) and the admissibility of 
the intrinsic classification algorithm (Theorems [5] and [6]). 

• In Section [21 we prove Theorem [3l in particular, we develop all the repre- 
sentation theoretic machinery that is needed for the proof. 

• Finally, in Appendix |Al we give the proofs of all technical statements which 
appear in the previous sections. 

Acknowledgement: The first author would like to thank Joseph Bernstein for 
many helpful discussions concerning the mathematical aspects of this work. He 
also thanks Richard Askey for his valuable advice about Legendre polynomials. The 
second author is partially supported by Award Number R01GM090200 from the Na- 
tional Institute of General Medical Sciences. The content is solely the responsibility 
of the authors and does not necessarily represent the official views of the National 
Institute of General Medical Sciences or the National Institutes of Health. This 
work is part of a project conducted jointly with Shamgar Gurevich, Yoel Shkolnisky 
and Fred Sigworth. 

1. Preliminaries and main results 

1.1. Setup. Let {V, (•, •)) be a three-dimensional, oriented, Euclidean vector space 
over R. The reader can take V = equipped with the standard orientation and 
(•, •) to be the standard inner product. Let W = CV denote the complexification 
of V. We equip W with the Hermitian product (•, •) : W x W ^ C, induced from 
(•, •), given by 

{u + iv, u + iv') ~ (m, v) + (f , v') — i (u, v') + i {v, u ) . 



THE CLASS AVERAGING PROBLEM 



9 



Let SO {V) denote the group of orthogonal transformations with respect to tlie 
inner product (•, •) which preserve the orientation. Let S {V) denote the unit sphere 
in V, that is, S {V) = {v eV : {v, v) = 1}. Let X = Fr {V) denote the manifold of 
oriented orthonormal frames in V , that is, a point a; G X is an orthonormal basis 
X = (ei, 62, 63) of V compatible with the orientation. 

We consider two commuting group actions on the frame manifold: a left action 
of the group SO iV), given by 

g t> (61,62,63) = (561,362,^63) , 

and a right action of the special orthogonal group 50(3), given by 

(61,62,63) < 5 = (01161+02162 + 03163, 

01261 + 02262 + 03263, 

01361 + 02362 + 03363), 

for 

(Oil 012 Oi3 
O2I 022 O23 
031 032 O33 

We distinguish the copy of 50(2) inside 50(3) consisting of matrices of the form 

(On 012 
021 022 
1 

and consider X as a principal 50 (2) bundle over 5 (V) where the fibration map 
■K : X ^ S {V) is given by tt (61,62,63) = 63. We call the vector 63 the viewing 
direction. 

1.2. The Transport data. Given a point v G S (V), we denote by Xy the fiber of 
the frame manifold laying over v, that is, Xy = {x £ X : n (x) = v} . For every pair 
of frames x,y G X such that tt (x) ^ ±7r (y), we define a matrix T {x, y) € SO (2), 
characterized by the property 

x<T{x,y) = t^{x),n{y) iv) ' 
where t7r(rE).7r(j/) • -^v^y) ^ Xtt{x) IS the morphism between the corresponding fibers, 
given by the parallel transport mapping along the unique geodesic in the sphere 
5 (V) connecting the points tt (y) with n {x). We identify M? with C and consider 
T(x,y) as a complex number of unit norm. The collection of matrices {T{x,y)} 
satisfy the following properties: 

• Symmetry: For every x,y G X, we have T (x, y) = T {x, y)~^, where the 
left hand side of the equality coincides with the complex conjugate T {x, y). 
This property follows from the fact that the parallel transport mapping 
satisfies: 

*7r(2/),7r(x) = *,r(a:),7r(y)- 

• Invariance: For every x,y £ X and element g G 50(F), we have that 
T [g <> x,g t> y) = T{x,y). This property follows from the fact that the 
parallel transport mapping satisfies: 

tiv{g>x),Tr{gt>y) = 9 ° in{x),ir{y) ° 9 ^ i 

for every g G SO {V). 
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• Equivariance: For every x,y £ X and elements gi, 52 G SO (2), we have 
that T {x <i gi,y <i §2) — gi^T (x, y) 52- This property follows from the fact 
that the parallel transport mapping satisfies: 

t7T(x<gi),7r{y<g2) ~ ^■k(x) ,'ir{y) ^ 

for every gi,g2 & SO{2). 
The collection {T(a;, y)} is referred to as the transport data. 

1.3. The parallel transport operator . Let H—C{X) denote the Hilbertian 
space of smooth complex valued functions on X (here, the word Hilbertian means 
that H is not complete jl, where the Hermitian product is the standard one, given 

by 



(/l, /2)« ^ J (^) /2 {x)dx, 

for every /i,/2 S 'H, where dx denotes the normalized Haar measure on X. In 
addition, H supports a unitary representation of the group SO {V) x SO (2), where 
the action of an element g = (^1,(72) sends a function s e "H to a function g • s, 
given by 

(g ■ s) (x) = s {g^^ [> X <J g2) , 

for every x £ X. 

Using the transport data, we define an integral operator T :% ^ % as 



T{s) (x) = J T{x,y)s{y)dy, 

for every s € H. The properties of the transport data imply the following properties 
of the operator T: 

• The symmetry property implies that T is self adjoint. 

• The invariance property implies that T commutes with the SO (V) action, 
namely T {g ■ s) = g ■ T (s) for every s e "H and g S SO {V). 

• The implication of the equivariance property will be discussed later when 
we study the kernel of T. 

The operator T is referred to as the parallel transport operator. 

1.3.1. Localized parallel transport operator. The operator which arise naturally in 
our context is a localized version of the transport operator. Let us fix an angle 
a £ [0, tt] , designating an opening angle of a spherical cap on the sphere and consider 
the parameter h = 1 — cos(a), taking values in the interval [0, 2]. 

Given a choice of this parameter, we define an integral operator Th :% % as 



(1.1) Th{s){x)^ J T{x,y)s{y)dy. 

y(EB{x,a) 

where B {x,a) — {y E X : (tt (x) , tt (y)) > cos (a)}. Similar considerations as before 
show that Th is self-adjoint and, in addition, commutes with the 5*0 {V) action. 
Finally, note that the operator Th should be considered as a localization of the 
operator of parallel transport discussed in the previous paragraph, in the sense, that 



'^In general, in this paper, we will not distinguish between an Hilbertian vector space and its 
completion and the correct choice between the two will be clear from the context. 




Figure 4. The first four eigen values of the operator T^, presented 
as functions of /i G [0, 2]. 



now only frames with close viewing directions interact through the integral (jl.ip . 
For this reason, the operator Th is referred to as the localized parallel transport 
operator. 

1.4. Spectral properties of the localized parallel transport operator. We 

focus our attention on the spectral properties of the operator Th, in the regime 
/i <C 1, since this is the relevant regime for the class averaging application. 

Theorem 3. The operator Th has a discrete spectrum A„(ft.), n G N, such that 
dimH (A„ (/i)) = 2n + 1, for every h G (0,2], Moreover, in the regime h <^ 1, the 
eigenvalue A„ (h) has the asymptotic expansion 

For a proof, see Section [2j 

In fact, each eigenvalue A„ (h), as a function of h, is a polynomial of degree n + 1. 
In Section [21 we give a complete description of these polynomials by means of a 
generating function. To get some feeling for the formulas that arise, we list below 
the first four eigenvalues 
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The graphs of Ai (/i), i ~ 1, 2, 3, 4 are given in Figure SI 
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1.4.1. Spectral gap. Noting that A2 (h) attains its maximum at ft. = 1/2, we have 

Theorem 4. For every value of h d [^i'Aj maximal eigenvalue of Th is Xi (h). 
Moreover, for every value of h d [0, 1/2], there is a spectral gap G{h) of the form 

G(h) = Ai {h) - A2 {h) = ift2 _ 1 ;j3^ 

2 D 

For a proof, see Appendix El Note that the main difficulty in proving the second 
statement is to show that A„ (/i) < A2 {h) for every h e [0, 1/2], which looks evident 
from looking at Figure |4l 

Consequently, in the regime /i ^ 1, the spectral gap behaves like 

G{h)-\h\ 

1.5. Main algebraic structure. We proceed to describe an intrinsic model W of 
the Hermitian vector space W , that can be computed as the eigenspace associated 
with the maximal eigenvalue of the localized parallel transport operator T^, pro- 
vided ft. ^ 1. Using this model, the Euclidean inner products between the viewing 
directions of every pair of orthonormal frames can be computed. 

• Extrinsic model: for every point x £ X, let us denote hy W 
be the unique complex morphism sending 1 G C to the complex vector 

/i - if2 e W. 

• Intrinsic model: we define W to be the eigenspace of Th associated with 
the maximal eigenvalue, which by Theorems [3] and |4l is three-dimensional. 
For every point x £ X , there is a map 

= ^2/3 • (eu^rlW)* : C ^ W, 

where ev^ : H — > C is the evaluation morphism at the point x, namely, 

eWa; (/) = ./ [x] , 

for every f G T-l. The pair (W, {cp^ : x £ X}) is referred to as the intrinsic 
model of the vector space W. 

The algebraic structure that underlies the intrinsic classification algorithm is the 
canonical morphism 

defined by 

t{v) {x) = Vm-Sl (v), 

for every x £ X. The morphism r induces an isomorphism of Hermitian vector 
spaces between W equipped with the collection of natural maps {S^ : C W} and 
W equipped with the collection of maps {ip^ : C — >■ W}. This is summarized in the 
following theorem: 

Theorem 5. The morphism t maps W isomorphically, as an Hermitian vector 
space, onto the suhspace W C "H. Moreover, 

for every x £ X . 
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For a proof, see Appendix |A] (the proof uses the results and terminology of 
Section [2]). 

Using Theorem [Sj we can express in intrinsic terms the inner product between 
the viewing directions associated with every ordered pair of frames. The precise 
statement is 

Theorem 6. For every pair of frames x,y Cz X, we have 

(1-2) {^{x),^{y)) = \{^Av).^y{n))\'l. 

for any choice of complex numbers w, u G C of unit norm. 

For a proof, see Appendix [Xj Note that substituting w = m = 1 in (|1.2I) we 
obtain dofg)) . 

1.6. Explanation of Theorems [l] and [2l We end this section with an explana- 
tion of the two main statements that appeared in the introduction. The explanation 
is based on inspecting the limit when the number of images N goes to infinity. Pro- 
vided that the corresponding frames are independently drawn from the normalized 
Haar measure on X (empirical assumption); in the limit: the transport matrix Tm 
approaches the localized parallel transport operator Th : V. ^ H, for some small 
value of the parameter h. This implies that the spectral properties of T/v for large 
values of N are governed by the spectral properties of the operator T/,, when h lies 
in the regime h ^1. In particular, 

• The statement of Theorem [1] is explained by the fact that the maximal 
eigenvalue of T\ has multiplicity three (see Theorem[3l) and that there exists 
a spectral gap G{h) ~ h/2, separating it from the rest of the spectrum (see 
Theorem |4]). The later property ensures that the numerical computation 
of this eigenspace makes sense. 

• The statement of Theorem [5] is explained by the fact that the vector space 
Wtv is a numerical approximation of the theoretical vector space W and 
Theorem [5] 

2. Spectral analysis of the localized parallel transport operator 

In this section we study the spectral properties of the localized parallel transport 
operator T/j, mainly focusing on the regime /i <C 1. But, first we need to introduce 
some preliminaries from representation theory. 

2.1. Isotypic decompositions. The Hilbert space H, as a unitary representation 
of the group 5*0(2), admits an isotypic decomposition 

(2.1) H=^Hk, 

where a function s £ Hk if and only ii s {x <J g) — g^s (a;), for every x £ X and 
g G 5*0(2). In turns, each Hilbert space H^, as a representation of the group 
SO{V), admits an isotypic decomposition 

(2.2) ^fc= 

nGN>0 

where 'Hn,k denotes the component which is a direct sum of copies of the unique 
irreducible representation of SO{V) which is of dimension 2n + 1. A particularly 
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important property is that each irreducible representation which appears in ()2.2p 
comes up with muhiphcity one. This is summarized in the foUowing theorem: 

Theorem 7 (Muhiphcity one). If n < \k\ then Hn,k = 0. Otherwise, Hn,k is 
isomorphic to the unique irreducible representation of SO{V) of dimension 2n + l. 

for a proof, see Appendix [X] 

The fohowing proposition is a direct imphcation of the equivariance property 
of the operator and foUows from Schur's orthogonahty relations on the group 
50(2): 



Proposition 1. We have 



CkerT;,. 

k^-l 

Consequently, from now on, we will consider as an operator from H-i to 
H^i. Moreover, since for every n > 1, Hn,-! is an irreducible representation of 
SO{V) and since Th commutes with the group action, by Schur's Lemma Th acts 
on Hn,-! as a scalar operator, namely 

n\nn^^l^Xn{h) Id. 

The reminder of this section is devoted to the computation of the eigenvalues 
A„ {h). The strategy of the computation is to choose a point xq & X and a "good" 
vector Un & T-Ln,-i such that m„ (a;o) ^ and then to use the relation 

Th {Un) [xq) = Xn (h) Un (xq) , 

which implies that 

(2.3) A„(M = ^iiMM. 

Un (Xo) 

2.2. Set-up. Fix a frame xq G X, xq — (ei, 62, 63). Under this choice, we can safely 
identify the group SO {V) with the group 5*0(3) by sending an element g S SO{V) 
to the unique element h g 50(3) such that g \> xq ~ xq <l h. Hence, from now 
on, we will consider the frame manifold equipped with a commuting left and right 
actions of SO (3). 

Consider the following elements in the Lie algebra so (3) : 

Ai = 



A. = 




A, = 



The elements A,;, i = 1, 2, 3 satisfy the relations 

[^3,^1] = A2, 
[^3,^2] = -Ai, 
[Ai,A2] = A3. 
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Let {H, E, F) be the foUowmg sh triple in the complexified Lie algebra Cso (3): 

H = -~2iA3, 
E = iA2-Ai, 
F = Ai+iA2. 

Finally, let {H^,E^,F^) and [H^.E'^.F^) be the associated (complexified) 
vector fields on X induced from the left and right action of 5*0(3) respectively. 

2.2.1. Spherical coordinates. We consider the spherical coordinates of the frame 
manifold w : (0, 2tt) x (0, tt) x (0, 27r) — > X, given by 

io{ip,e,a) =xo<\ e'^^^e^^^e"^^ 

We have the following formulas 

• The normalized Haar measure on X is given by the density 

sin {d) , ,^ , 
TrdiodOda. 

2{2Tif 

• The vector fields [H^,E^,F^) are given by 

= 2id^, 

E^ = -e-"f{ide + coi{e)d^-\/sm{e)do:), 
F^ = -e''P{ide-cot{e)d^ + llsm{e)dc,). 

• The vector fields (^H^ , E^, F^^ are given by 

= ~2ida, 

E" = e'" {ide + cot {9) d^-Xj sin (d) d^) , 
F^ = e''" {ide + cot {9) 9„ - 1/ sin [9) d^) . 

2.3. Choosing a good vector. 

2.3.1. Spherical functions. Consider the subgroup T C SO (3) generated by the 
infinitesimal element A^. For every k (z "L and n > k, the Hilbert space 'Hn,k 
admits an isotypic decomposition with respect to the left action of T: 

n 

m— — n 

where a function s S H™^. if and only if s (e~*"^^ \> x) = e*™*s (x), for every x € X. 
Functions in W^k usually referred to in the literature as (generalized) spherical 
junctions. Our plan is to choose for every n > 1, a spherical function w„ G _i 
and exhibit a closed formula for the generating function 

n>l 

Then, we will use this explicit generating function to compute w„ (a;o) a-nd Th (wn) (a;o) 
and use ()2.3p to compute A„ [h). 
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2.3.2. Generating function. For every n > 0, let ip„ G 'H'!^ g be the unique spherical 
function such that ip^ {xq) 1. These functions are the well known spherical 
harmonics on the sphere. Define the generating function 

n>0 

The following theorem is taken from [5]: 
Theorem 8. The function Go,o admits the following formula: 

Go,o iv^, 0, a, t) = {l- 2t cos {6) + t^y^^^ . 

Take m„ — E'^F^^'n- Note that indeed u„ e _i and define the generating 
function 



n>l 



It follows that, Gi__i = F^Gq q. Direct calculation, using the formula in 
Theorem [51 reveals that 

(2.4) Gi,-i{ip,e,a,t) = e~*("+'^)[3sin(0)'i2(i_2icos(0) + i')"'^/' 

-t cos (6*) (l - 2t cos (9) + t^) '^'^ 
-i(l-2icos {e) + t^ 



It is enough to consider Gi,_i when = a = 0. We use the notation Gi__i (0, t) 
Gi^_i(O,0,O,i). By (EH) 



(2.5) Gi,_i(6',t) = 3sin(6i)^i2^1-2icos(6l) + <2)"^/- 

-t cos (61) (1 - 2i cos {6) + t^y^'^ 

-^(l-2^cos(0) + ^2)-^/^ 

2.4. Computation of u„ (xq). Observe that 

Gi,_i(0,t) = Y,Un {xo)t^. 

n>l 

Direct calculation reveals that 

Gi,-i(0,i) = -2til-ty^ 

n>l ^ ^ 

Since {^^-i) ^ '^"'2 ("^ -'^)' obtain 

(2.6) u„ (xo) = -n (n + 1) . 
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2.5. Computation of Th (un) (a;o)- Recall that h — I — cos (a). 
Using the definition of Th , we obtain 

Th{un){xo)= J T {xo,y)u,i{y)dy. 

y£B{xo,a) 

Using the spherical coordinates, the integral on the right hand side can be written 

as 

1 f , f sin (9) 



{2nf 



dip J — - — d0 J T {xo,u} {ip,6,a))uniu} {ip,9,a)) . 





First 

(2.7) T{xo,uj{ip,9,a)) = T {xo,xq <J e'^^'e^'^^e"'^') 



where the third equality uses the invariance property of the transport data and the 
second equality uses the equivariance property of the transport data. 
Second, since u„ € "H^ _^ we have 







(2.8) u„ (u; {if, 9, a)) = e-'^^Un {xo < e'^-) e'^". 
Combining (12. 7p and (|2.8p . we conclude 

a 

(2.9) Th{un){xo) = I ^^^T{xo,xo^e''^^)un{xo<3e''^^)d9 

^«„(xo<e«-4^)d^. 



where the second equality uses the fact that xq < e^^^ is the parallel transport of 
xq along the unique geodesic connecting tt (xq) with tt (xq <\ e^"*^). 
Denote 

a 

In{h)^ f ^^u,,{xo<ie'^^)d9. 



2 

Define the generating function / {h, t) = J2 (h) and observe that 







n>0 





Direct calculation reveals that 
(2.10) I{h,t) = l/2[h{l + 2tih-l) + t^y^^^ 

-th (2 -h){l + 2t {h - 1) + t^y^^^ 
-t-\{l + 2t{h-l) + t'Y'^ - (1 - t))]. 
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2.6. Proof of Theorem \S[ Expanding I {h,t) with respect to the parameter t 
reveals that the function /„ (h) is a polynomial in h of degree n + 1. Then, using 
Equation (12.31) . we get 

n(n + 1) 

In principle, it is possible to obtain a closed formula for A„ (h) for every n > 1. 

2.6.1. Quadratic approximation. We want to compute the first three terms in the 
Taylor expansion of A„ (h): 

A„ (h) = A„ (0) + dhXn (0) + ^hh^ + o {h') . 



We have 



A„ (0) 
dhXn (0) 
dlX„ (0) 



4(0) 



Observe that 



n {n - 


f 1)' 


dhIn 


(0) 


n {n - 


f 1)' 




(0) 


n [n - 


f !)■ 


dlln 


(0). 



n>l 

Direct computation, using Formula (|2.10l) . reveals that 
/(0,i) = 0, 
9,J(0,t) = -;^n(n+l)r, 

dll{Q,t) - i5]r^(r^ + l)(l + (,^ + 2)(n-l))^^ 

n>l 

Combing all the above yields the desired formula 

A„(/.)4/.-i±i^i±|M^/.^ + 0(/.3). 
This concludes the proof of the theorem. 

Appendix A. Proofs 
A.l. Proof of Theorem [4l The proof is based on two technical lemmas. 

Lemma 1. The following estimates hold: 

(1) There exists hi e (0,2] such that A„ (/i) < Ai [h), for every n > 1 and 
h e [0,/ii]. 

(2) There exists /12 G (0,2] such that \n{h) < A2 (/i), for every n > 2 and 
h e [0,/i2]. 

The proof appears below. 

Lemma 2. The following estimates hold: 

(1) There exists Ni such that Xn (h) < Ai (h), for every n> Ni and /i G [/ii, 2]. 
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(2) There exists N2 such that A„ (h) < Ai (h), for every n > N2 and h G 

[h2,im. 

The proof appears below. 

Granting the vahdity of these two lemmas we can finish the proof of the theorem. 

First we prove that A„ (h) < Ai (h), for every n > 1 and h G [0, 2]: By Lemmas 
ll|2l we get that A„ (h) < Ai (h) for every h € [0, 2] when n > Ni. Then, we verify 
directly that An (h) < Ai (h) for every h e [0, 2] in the finitely many cases when 
n < Ni. 

Similarly, we prove that A„ (h) < A2 {h) for every n > 2 and h E [0, 1/2]: By 
Lemmas I1I2I we get that A„ (h) < A2 (h) for every h € [0, 1/2] when n > N2- Then, 
we verify directly that A„ (h) < Ai (h) for every h S [0, 1/2] in the finitely many 
cases when n < 

This concludes the proof of the theorem. 

A. 2. Proof of Lemma [TJ The strategy of the proof is to reduce the statement to 
known facts about Legendre polynomials. 

Recall h = \ — cos (a). Here, it will be convenient to consider the parameter 
z = cos (a) , taking values in the interval [— 1 , 1] . 

We recall that Legendre polynomials P„ (z), n G N appear as the coefficients of 
the generating function 



P{z,t) = {t^ - 2tz + l) 

Let 



-(1/2) 



2(T3I) " = 

(z) = { 2(T3IJ n = l . 
dz\n-i (z) n>2 
Consider the generating function 

00 

J(z,i) = ^ J„(z)i". 
The function J {z,t) admits the following closed formula 

/AiN T/ ,X t + tZ + t'^ + l,n ^, -,\-l/2 

(A.l) J{z,t)^ 2(l-z) +^^^ + ^) 

Using (jA.l[) . we get that for n > 2 

2(1 - z) 

where Qn (z) = (—1)" Pn (z). In order to prove the lemma, it is enough to show that 
there exists zq G (—1, 1] such that for every z £ [— l,zo] the following inequalities 
hold 

• Qn (z) < Qs (z) for every 71 > 3. 

• Qn (z) < Q2 (z) for every n > 2. 

• Qn (z) < Qi (z) for every n > 1. 

• Qn {z) < Qo (z) for every n > 0. 

These inequalities follow from the following technical proposition. 

Proposition 2. Let uq e N. There exists zq E (—1, 1] such that Qn (z) < Qn,, (z), 
for every z E [— 1, Zq] and n> ng. 
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The proof appears below. 

Take ho = hi = 1 + zq. Granting Proposition lA.ll verify that J„ (z) < J2 (z), 
for n > 2, z e [— l,zo] which imphes that A„ (h) < Ai (h), for n > 2, /i e [0, /iq] 
and Jn (z) < J3 (z), for n > 3, 2; € [—1, 20] which imphes that A„ (h) < A2 (h), for 
n > 3, /i e [0,ho]. 

This concludes the proof of the Lemma. 

A. 2.1. Proof of Proposition\^ Denote by ai < 02 < ... < a„ the zeroes of (5„ (cos (a)) 
and by /i]^ < < ... < the local extrema of Qn (cos (a)). 

The following properties of the polynomials Qn are implied from known facts 
about Legendre polynomials (Properties 1 and 2 can be verified directly), that can 
be found for example in the book PT: 

Property 1: < ji^ < a^+i, for i ~ 1, ..,n — 1. 

Property 2: g„ (-1) = 1 and d,Qn+i (-1) < d,Qn (-1) < 0, for neN. 
Property 3: \Qn (cos(^J)| > |(3„ (cos for i = 1, .., [n/2]. 

Property 4: (i - l/2)7r/r7, < a, < i7r/(n + 1), for i = 1,.., [n/2]. 
Property 5: sin(a)"^^^ • |(3„ (cos(a))| < ^Ij-Kn, for a e [0,7r]. 

Granting these facts, we can finish the proof. 
By Properties 1,4 

TT 27r 
2n '^^ n + 1 
We assume that n is large enough so that, for some small e > 

sin (ai) > (1 — e) fli. 

In particular, this is the situation when uq > N, for some fixed N — Ng. 
By Property 5 

|(3„ (cos(^i))| < yJlj-Kn ■ sin(^i)"^^^ 

< y/2/'Kn ■ sin (ai)^^^^ 

< ^/2/^-{{l-e)ai)-'/^ - 



Let flo G (0, vr) be such that Qng (cos (a)) > 2/Try/l — e, for every a < ag. Take 
zq = cos (ao). 

Finally, in the finitely many cases where no < n < N , the inequality Qn (z) < 
Qno (z) can be verified directly. 

This concludes the proof of the proposition. 

A. 3. Proof of Lemma [2l We have the following identity: 

(A.2) tr{T',) = '^, 

for every h e [0, 2]. The proof of (jA.2[) is by direct calculation: 

tr{Tl) = J T^ix,x)dx = 
xex 

fJ-Haar J Th (x, y) o Th (y, x) dx. 

x(£X yGB{x,a) 
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Since Th {x, y) = Th {y, x) ^ (symmetry property), we get 

a 

/^9\ f f , , f sm(d) ,„ 1 — COS (a) 
tr [Tl) = j j dxdy = j -^d6 = 

xeX yeB{x,a) 

Substituting, a = cos^^ (1 — h), we get the desired formula tr (T^) = h/2. 
On the other hand, 



(A.3) tr {Tl) =Y,tr (t^,^^ _ J = ^ + 1) ^" (^)' 



= 1 ' ■ n=l 

From ((X2|) and fO]) we obtain the following upper bound 

^/h 



(A.4) Xn{h)< 

vm + 2 

Now we can finish the proof. 

First estimate: We know that Ai (h) — h/2~h'^/8, hence, one can verify directly 
that there exists A^i such that Vh/y/An + 2 < Ai (h) for every n > Ni and h e 
[hi, 2], which implies by (IA.4[) that A„ (h) < Ai (h) for every n > Ni and h E [^i, 2]. 

Second estimate: We know that A2 (h) — h/2 — 5/i^/8 + h^ /6, therefore, one 
can verify directly that there exists N2 such that \/7i/y/4n + 2 < A2 (h) for every 
n > N2 and h E [h2,l/2], which implies by (|A.4p that A„ (h) < A2 (h) for every 
n > iV2 and /i e [/i2,l/2]. 

This concludes the proof of the Lemma. 

A. 4. Proof of Theorem [5l We begin by proving that t maps W — CV isomor- 
phically, as an Hermitian space, onto W = 71 (Amax {h))- 

The crucial observation is, that T-L (Amax (h)) coincide with the isotypic subspace 
V-i-i (see Section Consider the morphism a = ^2/3 ■ t : W ^ given by 

a{v) (x) =(5* (v). 

First claim is, that Ima C H-i, namely, that (5*^^ (v) — g^^S^. (w), for every 
V E W, X E X and g E 50(2). Denote by (•, ■)^^^ the standard Hermitian product 
on C. Now write 

{<>*x<giv) ,z)std = {v,S^<ig{z)) ^ {v,6^{gz)) 

= ((5* (v) , gz)std = {g^^S* (v) , z)std- 

Second claim is, that a is a morphism of 5*0 (V) representations, namely, that 
{9^) — ^g-^t>x ('^)j for every v E W, x E X and g E SO{V). This statement 
follows from 

(gv) , z) = {gv,Sa:{z)) = {v,g-^S^{z)) 

= {v, Sg-io-x (z)) = {S*g-i^^ (v) , Z),td- 

Consequently, the morphism a maps W isomorphically, as a unitary represen- 
tation of SO{V), onto Hi.-i, which is the unique copy of the three-dimensional 
representation of SO (V) in "H-i. In turns, this implies that, up to a scalar, a and, 
hence r, are isomorphisms of Hermitian spaces. In order to complete the proof it 
is enough to show that 

tr (r* o t) — 3. 
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This follows from 



3 

tr^TOT*) = —tr (a* o a) 

3 
2 

vGS{W) 

3 
2 

v£S{W) 

3 
2 

veS{W) xex 

3 



(a* o a (v) , w)^ dv 
{a (v) , a 

2dvdx = 3. 



2 

where dv denotes the normalized Haar measure on the five dimensional sphere 
SiW). 

Next, we prove that t o = ip^, for every x X. The starting point is the 
equation eVx\W o a = S*, which follows from the definition of the morphism a and 
the fact that Ima = W. This implies that o r = (5*. The statement now follows 
from 

^x°T = 0^ ip^O {t OT ) = d^OT 

This concludes the proof of the theorem. 

A. 5. Proof of Theorem [6l We use the following terminology: for every x G X, 
X = (ei, 62, 63), we denote hy Sx ■ C ^ V the map given by Sx {p + iq) = pei + qe2- 
We observe that Sx {v) = Sx (w) — iSx for every v £ C 

We proceed with the proof. Let x,y e X. Choose unit vectors Vx,Vy € C such 

that Sx (vx) — Sy (vy) = V. 
Write 

{Sx (vx) , Sy {vy)) = (Sx (vx) - iSx (ivx) ,Sy (vy) - iSy {ivy)^ 
(A.5) = {Sx {vx) ,Sy {vy)j + (^x {iVx) ,'5y (ivy)^ 

~i (^x (iVx) ,'Sy {vy)j + i (^x (vx) ,'Sy {ivy)j . 

For every frame z € X and vector e C, the following identity can be easily 
verified; 

Sz [ivz] = 7r(z) X Sz {vz) ■ 

This implies that 

SxiiVx) = IT (x) X Sx (Vx) ^ {x) X V, 

5y{ivy) = TT {y) X Sy{vy) ^ TT {y) X V. 
Combining these identities with Equation (jA.Sp . we obtain 

{Sx (vx) , Sy (vy)) = (v, v) + (tt (x) X v,tt {y) X v) - i (tt {x) xv,v)+i {v, tt {y) x v) . 
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Since v G ImSx H ImSy, it follows that {tt (x) x v,v) — {v,Tr{y) x v) — 0. In 
addition, 

(. (.) X . (,) X .) . dot l^^^^;!^^^^ ';)) = (.) , . (,)) . 

Thus, we obtain that {6x (vx) , Sy (vy)) = 1 + (tt (x) , tt (y)). Since the right hand 
side is always > it follows that 

(A.6) \{SAv.),Sy{vy))\ = l + inix),TT{y)). 

Now, notice that the left hand side of IA.6l does not depend on the choice of the 
unit vectors Vx and Vy. 

To finish the proof, we use the isomorphism r which satisfies t o 6x — ^x fo^' 
every x £ X, and get 

I (fx M,iPyivy))\^l + {T: (x) , TT (y)) . 
This concludes the proof of the theorem. 

A.6. Proof of Proposition O The basic observation is, that H, as a representa- 
tion of SO{V) X 5*0(3), admits the following isotypic decomposition 

oo 

where Vn is the unique irreducible representation of SO{V) of dimension 2n + 1, 
and, similarly, [/„ is the unique irreducible representation of SO(3) of dimension 
2n + 1. This assertion, principally, follows from the Peter Weyl Theorem for the 
regular representation of 5*0(3). 

This implies that the isotypic decomposition of Hk takes the following form 



where C/,^' is the weight k space with respect to the action 50(2) C 50(3). The 
statement now follows from the following standard fact about the weight decompo- 
sition: 

dim U„ = ' 



1 n > k ' 
This concludes the proof of the theorem. 
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